clear; clc; close all;

%% 固定参数
a0 = -2;
a1 = 0.83;
a2 = 0.14;
b0 = 1.3;
b1 = 0.1;
b2 = 1;
d1 = 0;
d2 = 0;

%% 不同初始 x0: [0, 4*pi, 15*pi, 19*pi]
x0_values = [-15*pi, 0, 4*pi, 19*pi];

% y0, z0 固定
y0 = 0.1;
z0 = 0.1;

%% 迭代设置
N_total = 10000;   % 总迭代步数
N_trans = 2000;    % 舍弃瞬态步数
N_show  = N_total - N_trans;

%% 准备颜色 (lines 函数返回一组可区分的颜色)
colors = lines(length(x0_values));

%% 创建 3D 图窗
figure;
hold on;
box on; grid on;
view(3);
xlabel('x_n','FontSize',12);
ylabel('y_n','FontSize',12);
zlabel('z_n','FontSize',12);
title('x_0 增强行为','FontSize',14);

for i = 1:length(x0_values)
    % 取当前初始 x0
    x_init = x0_values(i);
    col = colors(i,:);   % 对应一种颜色
    
    % 赋初值
    x = x_init;
    y = y0;
    z = z0;
    
    % 先做瞬态迭代
    for n = 1:N_trans
        x_new = x + y;
        y_new = sin(a0 * y^2 * sin(a1 + a2*x)) + b2*sin(z) + d1;
        z_new = b0*y + b1*z + d2;
        
        x = x_new; y = y_new; z = z_new;
    end
    
    % 记录后期轨迹
    X = zeros(1, N_show);
    Y = zeros(1, N_show);
    Z = zeros(1, N_show);
    
    for n = 1:N_show
        x_new = x + y;
        y_new = sin(a0 * y^2 * sin(a1 + a2*x)) + b2*sin(z) + d1;
        z_new = b0*y + b1*z + d2;
        
        x = x_new; y = y_new; z = z_new;
        
        X(n) = x;
        Y(n) = y;
        Z(n) = z;
    end
    
    % 绘制散点
    scatter3(X, Y, Z, 5, 'MarkerFaceColor', col, 'MarkerEdgeColor','none',...
        'DisplayName', sprintf('x_0=%.1f', x_init));
end

legend('Location','best');
